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Abstract 

Numerical simulations of active galactic nuclei (AGN) feedback in cool-core galaxy clusters have 
successfully avoided classical cooling flows, but often produce too much cold gas. We perform adaptive 
mesh simulations that include momentum-driven AGN feedback, self-gravity, star formation and stellar 
feedback, focusing on the interplay between cooling, AGN heating and star formation in an isolated 
cool-core cluster. Cold clumps triggered by AGN jets and turbulence form filamentary structures 
tens of kpc long. This cold gas feeds both star formation and the supermassive black hole (SMBH), 
triggering an AGN outburst that increases the entropy of the ICM and reduces its cooling rate. Within 
1-2 Gyr, star formation completely consumes the cold gas, leading to a brief shutoff of the AGN. The 
ICM quickly cools and redevelops multiphase gas, followed by another cycle of star formation/AGN 
outburst. Within 6.5 Gyr, we observe three such cycles. There is good agreement between our 
simulated cluster and the observations of cool-core clusters. ICM cooling is dynamically balanced 
by AGN heating, and a cool-core appearance is preserved. The minimum cooling time to free-fall 
time ratio typically varies between a few and > 20. The star formation rate (SFR) covers a wide 
range, from 0 to a few hundred M 0 yr -1 , with an average of ~ 40 M 0 yr -1 . The instantaneous SMBH 
accretion rate shows large variations on short timescales, but the average value correlates well with 
the SFR. Simulations without stellar feedback or self-gravity produce qualitatively similar results, but 
a lower SMBH feedback efficiency (0.1% compared to 1%) results in too many stars. 

Subject headings: galaxies: clusters: general, galaxies: intracluster medium, hydrodynamics 


1. INTRODUCTION 

X-ray observations of the intracluster medium (ICM) 
indicate that in the cores of some galaxy clusters, the 
gas with high density and low entropy has a cooling 
time much shor ter than the Hubble time (e.g. Fabian 
& Nulsen 1977). These clusters are often called “cool- 
core clusters. 77 In the absence of heating, cold gas is 
expected to condense out of the ICM in the center of 
cool-core clusters, driving an inflow (the so-called “cool¬ 
ing flow”) at rates of 100s-1000 M 0 yr -1 (see review by 
Fabian 1994]). However, more recent X-ray observations 
made by Chandra and XMM-Newton have revealed a 
dearth of cooler X-ray gas below 1 /3 of the virial temper ¬ 
ature (e.g. |Peterson et al. 2003 Peterson fc Fabian|2 006 ), 
and the observed star formation rate is usually much 


lower than the classical cooling ra te (e.g. McNamara & 
Q’Connell|1989 O’Dea et al.|2 008). The lack of observa¬ 


tional support for the existence of a classical cooling flow 
is referred to as the “cooling flow problem.” The solution 


ism that can balance radiative cooling 

(e.g., Ruszkowski 

& Begelman 2002f Ruszkowski et al.|| 

2004; Kim et al. 

2005; 

Heinz et al.|2006; Brighenti & Mathews||2006; Guo 

& Oh 

2008 Conroy & Ostriker 2008)). 
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is the major heating so urce in cool-core clusters (Mc- 


outbursts are often observed in nearby cool-core clusters 

( e -K-i 

Fabian et al. 

2006 

Wise et al. 

2007 

Baldi et al. 

2009; 

Blanton et a. 

201 

1), and the 

energy associated 


with t He cavities are usually enough to o ffset cooling 
(Dunn & Fabian 20061 Rafferty et al. 2006). There is 
mounting observational evidence linking the AGN activ¬ 
ities, the cooling properties of the ICM, and the mul¬ 
tiphase gas: clusters with low entropy X-ray gas (i.e., 
with a short cooling time) in the cores almost always har¬ 
bor line-emitting multiphase gas and radio-loud AGN, 


while the clusters with hotter cores usually do not (Best 


et al. 

2007; 

Cavagnolo et al. 2008; 

Birzan et al. 

IN umerical 

2012 

simuJ 

ations have shown 


in the form of momentum-driven jets can successfully 
suppress cooling in both cosmological sim ulations (e.g. 


Dubois et al. 

2010 Martizzi et 

al. 2012 

) and isolated 

systems (e.g. 

Gattaneo & Teyssier 2007| 

Gaspari et al. 

_ 1.11 ;_ L ' L 


the AGN is powered by the accretion of cold gas (Pizzo- 


lato fc Soker|2005 Gaspari et al. 2013| ) have not only"suc- 


cessfully achieved long-term global thermal equilibrium 
with a range of feedback efficiencies, but also produced 
multiphase filamentary structures (reaching radii of up 
to tens of kpc) that morphologically resemble the obser¬ 


vations of nearby cool-cor e clusters (e.g. |Gaspari et al. 
2012 Li fc Bryan| 2014a| hereafter LB14a). However. 


multiphase gas usually only exists at early times of the 
simulations (up to 1-2 Gyr), and due to the conserva¬ 
tion of angular momentum, much of the gas that has 
initially cooled forms a massive cold disk of 10 11 M 0 
that persists for many Gyr through the end of the sim- 
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illation (Li & Bryan 2014b hereafter LB 14b). Obser¬ 


vations oFoptical ana IK line emission shows great di¬ 
versity in the properties of multiphase gas in cool-core 
clusters: some systems host extended filamentary cold 


gas like the Perseus Cluster (Fabian et al. 2003, 2006) 


some only have detectio ns in the nu clei; some do not seem 


to harbor any cold gas ( Edwards et al.|[2009 McDonald 


et al.pOlO 


2011a). The amount of cold gas detected is 
also typically sm aller than 10 11 M 0 (Edge 2001 Salome 
& Combes 2003). This suggests that the phase with a 
long-lived massive cold disk is unrealistic and additional 
physics is needed in the simulations. 

Star formation is an important process that is often 
overlooked in simulations focused on suppressing cooling 
flows with AGN feedback. This is mainly because stellar 
feedback is energetically insufficient to bala nce cool ing 


except under extreme assumptions (Bregman & David 


19891 Skory et al. 

2013b. 

Sijacki &; Springe 

2006 


(e.g., 


star formation is suppressed but still occurs when AGN 
feedback is included. Numerous observations show that 
a significant fraction of brightest cluster galaxies (BC Gs^ 


19871 Cardiel et al. 

1995, 

Hicks & Mushotzkyf 

2005j ] 


dwards et al. |2007|. There 
is a general consensus that "a “reduced cooling flow : 


is 


present in most cool-core clusters: some of the hot ICM 
cools into (in some systems filamentary) cold structures; 
some of this cold gas forms stars, at a rate that is on av¬ 
erage an order of magnitude lower than the classical rate, 
but with a large variation, ranging from below detection 


per year (Hicks et al. 2010; 

McDonald et al. 2011b Hoffer 

et al. 

2012 ; McDonald et aJ 

. 2013). Star formation is also 

largel 

y responsible for ionizing the line-emitting filaments 

(O’Dea et al.|2010 McDonald et al. 2012). However, ex- 


the cooling of the ICM is not well understood, and there 
is a lack of detailed modeling of star formation in these 
largest galaxies. 

In this paper, we carry out three-dimensional adaptive 
mesh refinement (AMR) simulations to study the inter¬ 
play between ICM cooling, star formation and AGN feed¬ 
back in an idealized cool-core galaxy cluster. The model 
is based on our previous work that focused on the balance 
between cooling and momentum-driven AGN feedback 
(LB 14a & LB 14b). Besides cooling and AGN feedback, 
we also include the self-gravity of the gas, star formation, 
and stellar feedback in the standard run presented here. 
The goal of this work is to gain a better understand¬ 
ing of the following issues: (1) How does AGN suppress 
ICM cooling and mediate star formation? (2) What is 
the impact of star formation and stellar feedback on the 
AGN activities? (3) Given the self-regulating nature of 
AGN feedback, can we still put constraints on the feed¬ 
back efficiency? (4) What is the physical explanation for 
the observed relation between the ICM properties, AGN 
activities, the multiphase gas, and star formation? 

We structure the paper as follows: in Section [2j we 
describe the methodology of this work including the re¬ 
finement criteria, the setup of the initial conditions, and 
the modeling of AGN feedback, star formation, and stel¬ 
lar feedback; in Section [3j we describe the results from 


the standard simulation, analyze the intertwined rela¬ 
tionship between cooling, AGN and star formation, dis¬ 
cuss our parameter studies and other test runs we have 
preformed to understand the role of individual pieces of 
physics included in the simulations; in Section [4] we com¬ 
pare our results to the observations and other related 
simulation work, and discuss the potential effect of other 
physical processes such as thermal conduction. We con¬ 
clude this work in Section [5] 

2. METHODOLOGY 

Our three-dimensio nal simulations are performed using 
the AMR code Enzo (B ryan et al. 120141) , with the ZEUS 

: Norma 


hydrodynamic method ( [stone 


1992]). All of 


_ Borma n 

the simulations discussed in this paper have 64 3 root 
grids (N root ) and 10 maximum refinement level ( l m ax ) 
in a box of L = 16 Mpc. This gives the smallest cell 
size of Ax m i n = L/(N root 2 lmax ) ~ 244 pc. The resolu¬ 
tion is lower than LB 14a and LB 14b (with Ax m i n ~ 15 
and ~ 60 pc), but as is shown in LB 14a, the general 
results are converged as long as Ar m i n < 500 pc. The 
refinement strategy is kept the same as our previous work 
and is briefly o utline d here (a detailed description can be 
found in Li & Bryan 2012). A cell is refined if any of the 
following three criteria are met: (1) the gas mass in the 
cell exceeds 0.2 times the mean mass of the gas in one 
cell of the root grid; (2) the ratio of the cooling time to 
the sound-crossing time is smaller than some limit /?, (we 
use [3 = 6 here, a somewhat arbitrary value larger than 
1, to better resolve cooling); (3) the size of the cell is lar¬ 
ger than 1/4 of the Jeans length (the local Jeans length 
is always resolved by at least 4 cells in each spatial di¬ 
mension, to prevent nu merically -induced fragmentation 
of the gas) ( Truelove et al.|[l997 ). 

The important physical processes included in the sim¬ 
ulations are radiative cooling, momentum-driven AGN 
feedback and star formation (with feedback). The self¬ 
gravity of the gas is also taken into account. The radi- 


ative cool i ng cu rve is the same as that used in |Tasker 
& Bryan] (|2006|, with a temperature floor of 300 K. 
Throughout the paper, we refer to the diffuse ICM with 
temperatures above 10 7 K as the “hot” gas, and the con¬ 
densed gas with temperatures below 10 5 K as the “cold” 
gas. The simula tion data is analyzed using the yt pack¬ 
age (Turk et al.pOll). 

We describe the cluster setup in Section |2.1| The 
modeling o f AG N jets and star formation are described 
in Section |2.2|and Section |2.3[ respectively. 


2.1. Cluster Setup 


The initial setup is very similar to Li & Bryan (2012), 
LB 14a and LB 14b. Our idealized isolated cool-core 
galaxy cluster is constructed based o n the observation s 
of the Perseus Cluster. Following Mathews et al. (2006), 
the electron density profile is initially set to be: 


n e (r) = 


0.0192 


0.046 


0.0048 


l+(£t 


1 l.J 


1 + ( 200 ) 


1 1.1 


cm 


-3 


(1) 

where r is the distance to the cluster center in kpc. The 
power-law index of the last term is slightly steepened so 
that the density profile at large radii is more consistent 
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with cosmological simulatio ns as well as the observations 
of the outskirts of Perseus ([Urban et al.|[2014|). Because 
our focus is on the cluster core, and both cooling and 
dynamical timescales are very long outside the core, the 
exact value of the index is unimportant. 

The initial temperature prof ile with in r < 3 00 kpc is 
taken from observations (Ch urazov et al.| |2004): 


T = 7 


1 + (r/71) 3 
2.3 + (r/71)' 


keV 


( 2 ) 


while at larger radii, we a dopt the universal temperature 
profile found by Loken et ah; (2002), normalizing it to 
match the observations at r = 300 kpc: 


/ 3 r \ - 1 - 6 

T = 9.18(1 + --) keV, (3) 

V Zi T v i r / 

where r v i r = 2.440 Mpc is the virial radius of the cluster. 

The bri ghtest cluster galaxy ( BCG) is treated as a fixed 
potential (Mathews et al.||2006|) with the stellar acceler¬ 
ation: 


9*(r) = 


r 0.5975 

3.206 x lO " 7 


0.9 


r 1.849 

1.861 x 10- 6 


0.9' 


-1/0.9 


(4) 


The SMBH in the center of the cluster is treated as a 
point mass of Msmbh = 3.4 x 10 8 M© (Wilman et al.| 
2005). 

We assume an ideal gas law for the ICM with 7 = 5/3, 
and that the ICM is initially in hydrostatic equilibrium 
with the gravitational potential which includes the con¬ 
tribution from the dark matter, the BCG, the SMBH, 
an d the gas its elf. This allows us to fit an NFW profile 
(Navarro et al. 1996) to the dark matter halo: 


p(r) = 


Po 


±0 + ±Y 


(5) 


where po is found to be 8.42 x 10 -26 g cm -3 , and the scale 
radius R s = 351.7 kpc. 

We do not include any initial rotation or perturbation 
in the gas. 


2 . 2 . Jet Modeling 


The SMBH accretion and jet modeling is the same as 
in LB 14a & b. Here we only repeat the key aspects. 

We calculate the accretion rate Msmbh at each time 
step by dividing the total amount of cold gas in the close 
vicinity (r < 500 pc) of the SMBH by a typical accretion 
time (5 Myr). Then a fraction of the cold gas mass within 
the accretion region is removed in proportion to the cell 
mass. 

The AGN jets in our simulations are modeled as bi¬ 
polar outflows along t he z-axis from two parallel jet 


launching planes (Omma et al. 2004). The jet launch¬ 


ing planes are parallel to the x-y plane, at a distance of 
2 cell lengths from the SMBH located at the center of 
the simulation domain. Within the jet launching planes, 
the mass added to the cells at each time step At follows 


2 /r> 2 

Am oc e~ r 7 jey where r is the distance to the z -axis 
and Vj et = 1.5Ax m i n . Assuming only a small fraction of 
Msmbh is added onto the SMBH for radio mode AGN, 
the outflow rate is approximately equal to the accretion 


rate Msmbh- Therefore Am can be normalized by set¬ 
ting f Am = Msmbh At. Our simulations are far from 
being able to resolve the accretion physics. It is even 
possible that only a fraction of this accreted gas actually 
reaches the accretion disk, and the rest is pushed out by 
the jets. Strictly speaking, Msmbh is the rate at which 
gas is processed by the SMBH, but for simplicity, we refer 
to it as the SMBH accretion rate throughout the paper. 

The jet power is expressed as 


E = cMsmbhc 


(6) 


where e is the feedback efficiency. In our standard run 
we use e = 1%. We will discuss a test run with g 
lower feedback efficiency, e = 0.1%, in Section 3.3.4 
In all simulations performed here, we assume the jet is 
50% thermalized. Thus the kinetic power of the jet is 

^kinetic = E -^thermal = /kinetic^ = 0.5 E. As is shown 

in LB 14a, varying /kinetic from 0.1 to 1 does not signific¬ 
antly affect the results. To avoi d the energy loss throug h 
low-density channels (e.g., |Vernaleo & Reynolds||2006|, 
we enforce a small angle (0 = 0.15) precession on the jet 
with a period of 10 Myr. 


2.3. Star Formation and Stellar Feedback 
The star f ormation is modeled following |Cen fc Os- 


triker (11992 ). A star particle is created in a cell when 


the following set of criteria is met: 1 ) the gas density 


exceeds a certain value (1.67 x 10 - 24 gcm -3 is used in 
our simulations), 2 ) the gas mass exceeds the local Jeans 
mass, 3) the flow is convergent, and 4) the cooling time 
of the gas is shorter than the dynamical time of the gas 
in the cell (td yn = (37r/32 Gpceu) 1 ^ 2 ). When all of the 
criteria are met, a fraction of the gas within the cell is 
converted into a star particle: 


At o 

m * = e S F 7 - PgasAx 

^dyn 


(7) 


where csf is the star formation efficiency and is set to 
be 0.02 in our simulations. If m* is above the minimum 
mass for star particles, ra* ?m ^ n , then the star particle 
is created. When m* < a star particle is cre¬ 

ated with a probability of ra*/ra* 5m * n , and its mass 
would be m^^min or 80% of the gas mass in the cell, 
whichever is smaller. This minimum mass is set to 
prevent the formation of too many star particles which 
would slow down the computation, and our standard run 
uses m*, min = 1O 6 M 0 . 

The star particle is created over a period r following: 


^ f to t-t 0 , t~t 0 ^ 

mAt) = ra* / — 5 — exp(- )dt , ( 8 ) 

Jo r T 

where r is the maximum of the dynamical time td yn and 

10 Myr. 

We also include the stellar feedback from mostly Type 

11 SNe by injecting mass and energy back into the gas 
over a dynamical time. 25% of the mass of the created 
stars is returned to the gas phase as stellar winds and SNe 
ejecta. The feedback in the form of thermal energy from 
each star particle is esNe^*c 2 , where esNe is set to be 
10 -5 . Both mass and energy feedback are injected locally 
to the cell that contains the star particle. We discuss 
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how stellar feedback and the choice of star formation 


parameters affect the results in Section [373 
3. RESULTS 


In this section we p resent the main results of the sim¬ 
ulations. Section 13.11 describes t he c luster evolution in 
the standard run, and in Section |T2| we examine the in¬ 
tricate relationship between cool ing, star formation, and 
AGN feedback. In Section 13.31 we discuss other simu¬ 
lations that we have performed to study the role each 
physical process plays and the impact of changing some 
of the important simulation parameters. 


3.1. Cluster Evolution 

The cluster evolution in the first few hundred Myr is 
similar to that in LB 14a & b: about 300 Myr the 
initial minimum cooling time of the gas) after we start 
the simulation, runaway cooling happens first in the very 
center of the cluster. The cold gas is accreted onto the 
SMBH and turns on the AGN jets; the jets then per¬ 
turb the ICM in a non-linear fashion, dragging out low 
entropy gas with short cooling time to larger radii, caus¬ 
ing it to cool into clumps along the jet path. A de¬ 
tailed analysis of this cooling process can be found in 
LB 14b. These cold, dense clumps have a lower pressure 
than the surrounding ICM, which drives more hot gas 
to cool onto them, but as the clumps move through the 
hot ICM they also suffer from shock heating, ram pres¬ 
sure stripping and Kelvin-Helmholtz instabilities, which 
can reduce their masses or even destroy them completely. 
The fate of individual clumps is uncertain (as they are 
numerical difficult to model accurately), but at the be¬ 
ginning of the clump formation stage the growth gener¬ 
ally overwhelms destruction. 

The cold clumps fall to the center of the cluster roughly 
within a dynamical time. Some of the cold gas is accreted 
onto the SMBH, boosting its power, which increases the 
entropy of the ICM, halting clump formation. Some of 
the cold clumps also form stars, and stellar feedback 
heats up the cold gas locally, reducing its mass. The 
total amount of cold gas in the first 1.5 Gyr peaks at 
only 7 x 10 9 M Q and does not form a rotationally sup¬ 
ported disk (the first row of Figure]]]). At t ^ 1.5 Gyr, all 
the cold gas has vanished and the AGN shuts off. The 
cluster enters a state that is very similar to the initial 
condition except that the gas is now turbulent. 

The second cycle is similar to the first one. At the 
end of the first cycle, the shortest cooling time of the 
ICM is only a few hundred Myr. With the AGN heat¬ 
ing off, the ICM quickly cools again. Due to the tur¬ 
bulent motion of the gas, runaway cooling first happens 
^15 kpc from the cluster center (the second row of Fig¬ 
ure [TJ from a blob of low entropy gas at t = 1.62 Gyr. 
The cooled gas shortly gets accreted onto the SMBH, 
turning the AGN back on. The AGN jets trigger more 
hot gas to cool into filamentary structures, resulting in 
a burst of star formation and AGN feedback (the third 
row of Figure [l]). The total amount of cold gas peaks 
at 5 x 10 10 M 0 , about one order of magnitude higher 
than the previous cycle. Cooling out of the ICM at a few 
tens of kpc, the cold gas carries a non-zero total angu¬ 
lar momentum. Due to the larger amount of cold gas in 
this cycle (and higher total angular momentum), SMBH 
accretion and star formation cannot consume all of it 


quickly enough before the gas settles to a rotating disk 
in the center of the cluster at t^ 2.5 Gyr. The cluster 
enters a stage that appears very similar to the last stage 
of the simulations with AGN feedback only (LB 14a): the 
AGN outburst has elevated the entropy of the ICM in 
the cluster core; all the gas with short cooling time (~ 
a few hundred Myr) has already cooled, and thus clump 
formation has stopped; the ICM cools directly onto the 
rotationally supported cold disk; the SMBH continues 
to accrete from the inner part of the cold disk, keeping 
the AGN jets on, which balances the ICM cooling (the 
fourth row of Figure [l]). In the simulations with AGN 
only, due to the conservation of angular momentum, a 
disk of ~ 10 11 M 0 lasts through the end of the simu¬ 
lation without showing signs of shrinking (LB 14a). In 
the standard run shown here, however, star formation 
gradually consumes the cold disk. By t = 3 Gyr, the 
mass of the disk has to be reduced to only r\j 10 10 M 0 . 
At t ~ 4.3 Gyr, the disk is completely gone (the last row 
of Figure [Tb . This marks the end of the second cycle. 

The third cycle begins at t = 4.8 Gyr and ends around 
t = 6 Gyr. At the end of our simulation at roughly 
t ~ 6.5 Gyr, there is no cold gas, and both AGN and 
star formation are off, but the ICM is cooling and ready 
to enter yet another cycle. 

Figure [I] also shows that the star forming regions (the 
fourth column) are spatially correlated with the cold 
dense gas (first two columns). A comparison between 
the morphology of the cold gas in our previous higher 
resolution simulations r] (LB 14afc b) an d the high quality 
UV images of BCGs m CLASH (jPostman et al. 2012) 
is presented in Donahue et al. 2015 (submitted to ApJ). 
Interestingly, the cycle we see in our simulation, in which 
cold, filamentary star-forming gas settles into a cent¬ 
ral disk and then disappears, also ec hoes the hyp othes- 
ized dust-cloud evolution sequence in Lauer et al. (2005) 
based on Hubble observations of 77 early-type galaxies. 


3.2. The Interplay between Cooling, AGN Feedback and 
Star Formation 

Both AGN and star formation are fueled by the cold 
gas that cools out of the hot ICM. AGN feedback heats 
the ICM and regulates its cooling rate, while star forma¬ 
tion consumes the cold gas and thus affects the AGN duty 
cycle. In this section, we study the intricate relation¬ 
ship between the ICM cooling, AGN feedback and star 
formation. There are many ways of measuring the cool¬ 
ing properties of the hot ICM in the core of the cluster, 
such as temperature, central entropy, and classical cool¬ 
ing rate (which is often defined as the cooling rate within 
the cooling radiu s where the cooling time is a few Gyr; 
e.g., Mit tal et al.|2 009). It is also found that precipitation 


in cool-core clusters is tightly linked to the minimum ra¬ 
tio of the cooling time o ver the free-fall time t coo \/t^ (e.g. 


Voit & Donahue| 2014|). Since central entropy, classical 
cooling rate” and the minimum ratio convey very sim¬ 
ilar information (low central entropy translates to high 
classical cooling rate and low minimum ratio), we mainly 


5 The cold clumps are better resolved with higher resolution and 
the bipolar cold structure along the jet path is more prominent at 
the beginning of the AGN outburst. Because we are not focused 
on the cold structure in this work, a lower resolution is chosen to 
reduce computational cost. 
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Fig. 1.— The projected gas density, temperature (weighted by density), X-ray luminosity, young stellar density and stellar density in the 
central 80 X 80 X 80 kpc 3 box at t = 0.56, 1.62, 1.73, 2.62 and 4.56 Gyr. The projection is along the rr-axis and the AGN jets are along the 
z-axis (vertical). The young stellar density only includes star particles less than 200 Myr old, which can be observed in the UV. The stellar 
density only shows the stars formed in the simulation but not the stars that are already in the BCG at the beginning of the simulation 
(treated only as a fixed potential). An animation of the temperature evolution can be viewed here: https://vimeo.com/115825854 
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use the minimum ratio to characterize the cooling prop¬ 
erties of the ICM and we show its evolution with time 
in Figure |2j along with the minimum cooling time t coo i, 
SMBH accretion rate, SFR and the total amount of cold 
gas. Note that we use the shell averaged t coo \ as is done 
in observations, instead of the actual value measured in 
each cell. 
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Fig. 2. — The time evolution of the minimum t C oolAff ratio, the 
minimum cooling time t c 00 i, the SMBH accretion rate MsmbHj 
the SFR, and the total amoun t of cold gas in the system. Data 
is sampled every 10 Myr and the color indicates time (from red to 
purple). The black line in the third panel is the Msmbh averaged 
with a 200 Myr moving window. 


The minimum t coo i/te ratio approaches its local min¬ 
imum at the beginning of every cycle. When AGN feed¬ 
back is just turned on, this ratio still decreases for a very 
brief period of time. It takes roughly one dynamical time 
(a few tens of Myr at r ^ 10 — 30 kpc) for most of the 
cold clumps to reach the SMBH and trigger the AGN 
outburst. The instantaneous SMBH accretion rate shows 
very large variations, but when we take its average with a 
moving window of 200 Myr, the trend is clear (the black 
line in the third panel of Figure [2]). The minimum ratio 
starts to increase right before the^vGN outflow rate (also 
equal to the SMBH accretion rate in our model) reaches 
its peak. In every cycle, within 1 Gyr of the time that 


the AGN is turned on, all the gas with short t coo i of a few 
hundred Myr has either cooled or been heated to higher 
temperatures by the AGN jets. The minimum t coo i/te 
ratio is elevated. In cases where star formation cannot 
consume the cold gas quickly enough and a rotationally- 
supported disk forms (e.g., the second cycle), this ratio 
plateaus because the AGN self-regulates through the cold 
disk (see LB 14a for details). As soon as the cold gas van¬ 
ishes, which turns off AGN feedback, the ICM starts to 
cool. As one would expect, precipitation begins again 
after a period ^ t c00 i, marking the beginning of the next 
cycle. 

The evolution of the star formation rate (SFR) is very 
similar to that of the average Msmbh, while the history 
of the total amount of cold gas basically mirrors the SFR. 

Both t coo \ and the minimum t coo \/t^ ratio indicate the 
general cooling properties of the ICM, but do not directly 
translate to the mass deposition rate - i.e., the rate at 
which the ICM is cooling into the cold phase. Since an 
AMR code does not trace the history of fluid elements, 
we can only estimate this mass deposition rate M CO oiing- 
We can measure the changing rate of the total amount 
of cold gas M C oid which is related to M coo ii n g- 

M co id = M coo ii ng — Msmbh — SFR — M heate d , (9) 

where the last term includes many physical processes 
that can reduce the amount of cold gas but are not 
strictly measurable in the simulations: shock heating, 
ram pressure stripping (of the cold clumps as they move 
through the hot ICM), Kelvin-Helmholtz instability (at 
the interface between cold clumps and the ICM) and SN 
heating. The first few processes are only significant when 
a large amount of cold clumps are present. Besides star 
formation, the depletion of the cold gas is dominated by 
SN heating. 

In an ideal situation where the thermal energy from 
SN feedback heats the cold gas exactly to the ICM tem¬ 
perature, we can write 


^sn = esNeSFR x c 2 = — SNheated x -k b TicM, (10) 

lim p 2 

where k b is the Boltzman’s constant and esNe = 10 -5 
(Section 2.3). Taking Ticm = 3 keV and assuming a 
typical delay of 20 Myr between star formation and SN 
feedback, we have 


M SNheated (t) = 1.23 x SFR(t - 20 Myr) . (11) 

In reality, SN might overheat a cell to higher temper¬ 
atures than the ICM temperature, resulting in an en¬ 
ergy “waste.” More often, it under-heats the cell such 
that the temperature goes up to only rsj 10 6 K or less 
where the cooling time is very short, and the cell quickly 
radiates the energy away and remains cold. Therefore, 
the calculation above (assuming optimum heating) likely 
overestimates the actual MsNheated in the simulation and 
should be taken as an upper limit. 

Approximating M heate d as MsNheated, we can rewrite 
Equation [9] as 

M coo ii ng ~ M co id + Msmbh + SFR + MsNheated • (12) 
Figure [3] shows the estimated mass deposition rate 
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Fig. 3. — The estimated mass deposition rate M C ooiing (light green), the SFR (blue) and the SMBH accretion rate Msmbh (grey) as a 
function of time. The solid green line and the black line show M C ooling and Msmbh averaged over 200 Myr, respectively. The dashed red 
line (300 Mq yr -1 ) shows the classical cooling flow prediction for M coo ii ng for the Perseus cluster. 



^cooling as a function of time. Again, Mcooiing shows 
a similar trend to the SFR and Msmbh, but peaks about 
20 Myr earlier at the beginning of every cycle. As a ref¬ 
erence, the red dashed line indicates the classical cooling 
rate of Perseus, which is recovered in our previous pure 


cooling flow simulation without AGN feedback (Li & 


Bryan 2012 LB 14a). M coo ii ng is significantly suppressed 
compared with the classical cooling flow prediction. 

Over a period of 6.5 Gyr, about 6 x 10 11 M 0 of 
gas has cooled in the simulation, with roughly half (< 
3 x 10 11 M 0 ) going into forming stars and the other half 
heated by supernovae (Figure [4|. The amount of gas 
that has been processed by the SMBH is 2.8 x 10 10 M o . 

Figure [5] shows that the radial profiles of the ICM 
density, temperature and pressure exhibit variations with 
time, but (outside of about 10 kpc) they never deviate 
much from the initial conditions and a cool-core appear¬ 
ance is preserved throughout the simulation (in partic¬ 
ular there is a central core of gas which is cooler than 
the bulk of the cluster). The ICM cooling is balanced by 
AGN heating in a dynamical way: the jet power over¬ 
whelms cooling during the first half of every AGN cycle; 
it gradually declines until cooling start to dominate and 
eventually leads to precipitation which triggers the next 
cycle (top panel of Figure |6|. At the end of the simu¬ 
lation, roughly 2/5 of the total jet energy is deposited 
within 100 kpc, and 3/5 within 300 kpc (bottom panel 
of Figure [6]). These numbers are comparable to those 
found in the AGN-only simulations (LB 14a). 

We have neglected the heating from SNe in Figure [6| 
SN feedback heats the local cold gas effectively, but its 
impact on the hot ICM is minimal. Our SN feedback 
recipe corresponds to 10 51 erg thermal energy output for 
every 180M o of stars formed, or 


Esn = 1-76 x 10 41 (j 


SFR 


M©/yr' 


ergs 


-l 


(13) 



Fig. 4.— The integrated amount of gas that is processed by the 
SMBH (grey), the gas that has formed stars (blue), the cold gas 
that exists in the system (green), and the estimated amount of cold 
gas that has been heated by SN feedback (red). 


The average SFR in our simulation is 40 M 0 yr 1 , so the 
average heating rate from SNe is < Esn >= 7 x 10 42 
erg s _1 , three orders of magnitude lower than the AGN 
heating rate. 


3.3. Other Simulations 

In this section, we discuss the other simulations that 
we have performed in order to understand the ro le that 
individual physic al proc esses play (Section |3.3.1 AGN 
feedback; Section |3.3.2| stellar feedback; Section |3.3T‘ 
self-gravity), and to study how the choices of simulation 
parameters affect the results. 





























































Fig. 5.— Gas density, temperature and pressure profiles weighted bv the X-ray emissivity. Profiles are plotted every 500 Myr since the 
AGN is first turned on, following the same color scheme as Figure p] The black line shows the initial condition at t = 0. The profiles 
bounce around the initial condition and a cool-core appearance is maintained. 
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Fig. 6.— Top: the AGN jet power (~the heating rate, shown in 
pink) and the total ICM cooling rate in the central r < 100 kpc 
(blue) and r < 300 kpc (green). The red line is the AGN heating 
rate averaged over 200 Myr. Bottom: The accumulated amount 
of heating from the AGN and the total cooling loss in the central 
r < 100 kpc (blue) and r < 300 kpc (green) of the cluster. 


3.3.1. No AGN feedback 

As a sanity check, we perform a simulation with star 
formation and stellar feedback but without AGN feed¬ 
back. As expected, a classical cooling flow develops 
within a few hundred Myr and stars are forming at hun¬ 


dreds of M 0 yr _1 . There is noticeable heating from SNe 
in the very center of the cluster as compared to a run 
without star formation or AGN feedback. However, the 
heating rate (estimated using equation 13) is only a few 
times 10 43 erg s -1 , far below the ICM cooling rate. 

This confirms that although star formation and stellar 
feedback have a significant impact on the cooling-AGN 
feedback cycle, AGN feedback still plays the major role 
in heating the ICM. 

3.3.2. No Stellar Feedback 

Prior to our main simulation, we also performed runs 
that include star formation but no stellar feedback. The 
results are very similar to our standard run. This is 
simply because stellar feedback is only roughly as effect¬ 
ive as star formation itself in reducing the amount of cold 
gas in the system (Figure [3. The primary difference we 
observe is that without stellar feedback, the star forma¬ 
tion rat e is higher than in our s tandard run, as one would 
expect ( Tasker fc Bryan||2006 ). 


3.3.3. The Effect of Self-gravity 

We have also performed simulations to test the import¬ 
ance of including the self-gravity of the gas. The grav¬ 
itational potential is dominated by the dark matter on 
large scales and by the stellar component of the BCG in 
the innermost ^10 kpc, while t h e co ntribution from the 
ICM is negligible (Li & Bryan 2012). This is the main 
reason why previous”imulations focusing on the AGN 
feedback usually ign ore the self-gravity of the gas (e.g., 
Gaspari et al. 2011 LB 14a & b). However, self-gravity 
may assist the development of gravitational instabilities 
of the cold disk, which can help redistribute angular mo¬ 
mentum and increase the accretion rate onto the SMBH 


(Lodato 2007). 


Our first test run includes AGN feedback and self¬ 
gravity but no star formation. The cluster evolution 
within the first 1.5 Gyr is very similar to the standard 
run discussed in LB 14a without the self-gravity of the 
gas. The cold gas settles into a rotating disk by t = 2 
Gyr. Self-gravity triggers gravitational instabilities, and 
the disk shows structures that resemble the spiral arms 
of disk galaxies. The accretion rate is higher compared 
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with LB 14a because the spiral-mode instabilities trans¬ 
port angular momentum, and enhances accretion. Thus 
the core temperature is also more elevated. The amount 
of cold gas is reduced, but the ICM is still cooling. The 
SMBH accretion rate is enhanced due to self-gravity, but 
it is still not enough to consume all the cold gas in the 
cluster center by the end of the simulation (within a few 
Gyr). Therefore, the problem of a long-lived massive 
cold disk persists, and no clump formation is seen after 
t = 2 Gyr in the over-heated cluster core. This means 
that self-gravity may alleviate the problem but is not the 
key to solving it. 

The other test we have performed is a simulation with 
all the ingredients included in our standard run except 
self-gravity. We find that the absence of self-gravity does 
not significantly alter the results. The cluster still exper¬ 
iences AGN outburst cycles. The main difference is that 
without self-gravity, the SFR is less bursty than in our 
standard run, peaking only at rsj 150 M 0 yr h This is 
likely because self-gravity can affect the density distribu¬ 
tion of t he star-forming clouds and speeds up the initial 
collapse (Kritsuk et al. |2011| |Lee et al.||2014|) . Another 
difference is that without the assistance of self-gravity, 
the accretion of the cold disk onto the SMBH is less effi¬ 
cient. Thus, within 6.5 Gyr the cluster only experiences 
two major cycles rather than the three seen in our stand¬ 
ard run. The cold disk also does not completely vanish 
(but is significantly reduced in mass) at the end the first 
cycle. 


3.3.4. SMBH feedback efficiency 

Previous simulations have found that due to self¬ 
regulation, AGN feedback with a wide range of feedback 
effi ciencies (e = 0.1%—1%) produces similar results (Gas- 
20121 LB 14a). This changes when star forma- 
d. 


pari et al. 


tion is inc 


udeo 


Although star formation and SMBH feedback evolve 
hand-in-hand (Figures [2] and [3| as the ICM goes through 
cycles of cooling and Treating, they also compete with 
each other at any given time as they share the same fuel 
- cold dense gas. In fact, as discussed earlier, the amount 
of cold gas that feeds star formation is more than an order 
of magnitud e hi gher than that accreted onto the SMBH 
(see Section 3.2 or Figure 13). 

In LB 14a, we chose e = 0.1% over e = 1% as our 
preferred value (and thus used this value in the stand¬ 
ard run) mainly because a feedback efficiency of e = 1% 
slightly overheats the core at late times. In the standard 
simulation discussed here, we use e = 1% because with 
most of the cold gas going into forming stars, a low AGN 
feedback efficiency is no longer sufficient. 

Figure [7] shows the cluster evolution for the simulation 
with e = 0.1% and everything else kept the same as the 
standard run. AGN feedback still manages to prevent the 
development of a classical cooling flow, and the cluster 
still shows somewhat cyclical behaviors similar to the 
standard run. This is not surprising since the average 
A^smbh here is roughly 10 times that in the standard run, 
making up for the lower efficiency, manifesting again the 
self-regulating nature of AGN feedback. However, the 
total mass of the cold gas in the system is almost always 
above 1O 1O M 0 and the SFR is also too high, with an 
average of 120 M 0 yr -1 . Thus a low feedback efficiency 
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Fig. 7. — The time evolution of the minimum t coo \/tff ratio, the 
minimum cooling time £ coo i, the SMBH accretion rate Msmbh> 
the SFR, and the total amount of cold gas in the simulation with 
a lower feedback efficiency e = 0.1%. The data is sampled every 
10 Myr and the color scheme is the same as Figure [2] The low 
efficiency results in too much star formation and too much cold 
gas. 


of e = 0.1% can be ruled out. 


3.3.5. Star Formation Prescription 

There are many parameters relat ed t o the star form¬ 
ation routine described in Section 12.31 The values we 
have chosen for our simulations are typically used in 
simulations of star forming galaxies that produ ce resu lts 
roughly cons i stent with the observation s (e.g., Tasker & 
Bryan 2006 Hummels & Bryan |2012). A systematic 
study of varying the values for all those parameters is 
beyond the scope of this work. Here we briefly discuss 
the results from a test run with a lower minimum mass 
for star particles (ra* jm i n = 5 x 10 4 M 0 , 1/20 of the value 
used in our standard run). 

Because star particles are easier to form at the begin¬ 
ning of the cooling cycle, more cold gas goes into form¬ 
ing stars and less gets accreted onto the SMBH, which 
leads to less energetic AGN outbursts. As a result, less 
condensation is triggered by the jets, and the first cycle 
is shorter (lasting for less than 1 Gyr) and gentler than 
that in the standard run. After global turbulence is built 
up by t = 1 Gyr, the cluster behaves more similarly to 
our standard run. The range of SFR and the minimum 
Gooi/te ratio stay roughly the same, except that the tail 
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of the SFR distribution extends a bit furth er, as one 
might expect. As discussed in Section |3.2[ the direct 
effect of star formation and stellar feedback on the hot 
ICM is negligible. We do not expect the results to change 
qualitatively when varying star formation parameters be¬ 
cause the battle is ultimately between ICM cooling and 
AGN feedback. However, due to the complex relations 
between cooling, star formation and AGN feedback, the 
exact choice of those parameters could have an effect on 
the strength and duration of the long term (> Gyr) AGN 
cycles. 


4. DISCUSSION 

In this section, we compare our results with observa¬ 
tions and other simulations done previously (Section |4.1| ), 
and we discuss the possible impact of other phy sical pro¬ 
cesses on the simulation results (in Section |4~2| ). 

4.1. Comparison with Observations 

A number of recent simulations show that gas con¬ 
densation occurs when the minimum t coo \/ts ratio of the 
system drops below a certain threshold. This threshold 
is fou nd to be ~ 1 in simula tions with idealized heat¬ 
ing in Mc Court et al. (2012 |), and ~ 10 in a spherica l 
geometry ( jSharma et a 20 I2l |Singh fc Sha rma||2015lh 
which is confirmed in |Gaspari et al. (|2U12|) with more 
realistic AGN feedback! LB 14 b finds that the value of 
this critical ratio can vary between ^3 — 10 depending 
on the strength of the perturbation and the way AGN 
jets interact with the ICM (dredging up gas with shorter 
£ coo i to larger radii where % is longer). 

Figure [8] shows the distribution of the minimum 
tcooi/tff ratio for the times when cold gas exists in the 
system. Precipitation can happen when the minimum 
t c ool/tff ratio is as low as ^ 1 (only very briefly) when the 
AGN is off, and as high as ^ 20 at the peak of the AGN 
outburst when the ICM is most perturbed. Over-plotted 
in Figure [8] are the observations of cool-core clusters with 
Ha detection. The bulk of the two distributions are con¬ 
sistent with each other. The observed distribution shows 
a longer tail at the high end and a deficit at the very 
low end where the minimum t coo \/t^ ratio is < 3. This 
discrepancy may be due to mergers and/or cosmological 
infall, effects that are not included in our isolated simu¬ 
lations. 

Observationally, in precipitating systems the radial 
profile of t coo i is usually ab ove ~ 10 % and never be¬ 
low 5 tff (Figure 1 of Vo it et al. 2014). Our simulation 
generally agrees with this finding (Figure [9|. Note that 
we plot only the initial t ff for clarity, but at later times 
the stars formed in the simulation also add mass and 
thus slightly steepen t ff. Despite this effect, our system 
is still found to have relatively low t coo i (between 5 and 
10 tff) more often than the observations seem to indicate. 
This may be partially due to the difference in the way 
t coo i is measured in simulations and in observations. For 
example, the yellow line that shows a dip at r ~ 20 kpc 
corresponds to the off-centered cooling at the beginning 
of the second cycle. Observationally, this off-centered 
low entropy blob could be identified as the X-ray center 
and thus change the spher icall y-averaged t coo 

As is shown in Section |3.2[ star formation is closely 
linked to ICM cooling, which can be characterized with 
the cooling time, t coo \. Observations confirm that they 



Min(t coo i/tff) 

Fig. 8. — The normalized distribution of the minimum t C ool/% 
ratio at times when cold gas exists in our standard run (blue) com¬ 
pare d with the observations (red) for precipitating systems (data 
from |Voit fe Donahue(2014)). 



r (kpc) 


Fig. 9.— The cooling time %>ol profile sampled every 500 Myr 
since the AGN was first turned on. The color scheme is the same 
as previous figures (from red-early to purple-late). The solid black 
line is 10% and the grey shaded region is 5% to 20%. The dashed 
black li ne is the conductively balanced solution and is discussed in 
Section 


are indeed correlated: Ha luminosity Lu a tends to be 
higher in systems with lower central entropy or smaller 
minimum t coo \/t^ ratio; it is often undetected if the ICM 
is hotter than s ome t hresh old (e.g. Ko > 30 keV cm 2 or 


is notter tnan s ome tnresnoiu je.g. x, ouKev cm or 
Gooi/tff > 20) ([Cavagnolo ~al. 2008; McDonald et al. 


2010 Rawle e t al.|2012 ' Voit &; Donahue|2014 ). We show 


in Figure TU the relation between the minimum t coo \/t ff 
ratio and the estimated Lu a in our simulation. To obtain 
I/Ha from star forming regions, we adopt the commonly 
used linear relat ionsh i p: L^ a (erg s s -1 ) = 1.26 x 10 41 
SFR (M© yr -1 ) (Kennicutt] 1998). Besides photoioniza¬ 
tion, some Ha emission may"also result from collisional 
ionization and the central AGN. For simplicity, we only 
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Fig. 10. — Left: the observed Ha luminosity I/Ha vs. the minimum £ C ool/Iff ratio in|Voit & Donahue| |2014). The black squares are the 
systems with Ha: detection, and the pink squares are the ones with only upper limits. Kight: the estimated Lh q vs. the minimum t coo lAff 
ratio in the standard run. The data is sampled every 10 Myr and the color scheme is the same as Figure [2] (from red-early to purple-late). 
The circles correspond to the I/Ha associated with star forming regions and the vertical bars indicate the estimated maximum contribution 
from collisional ionization. The shaded grey region in both panels denotes the range that the observed precipitating systems cover. 


calculate the emission from collisional ionization through 
i nterpolati on of a t able generated with Cloudy 13.00 
([Fbrland et al.|2013|) assuming coronal equilibrium. This 
gives an upperTimit because we are double-counting the 
photoionized regions. The contribution from collisional 
ionization is shown as the vertical bars and is much lower 
than photoionization most of the time (smaller than the 
size of the circles). 

Our system generally occupies the same area on the 
I/Ha-min(t coo i/tff) plot as the observed co ol-core clusters 


with detected Ha emission (data from Voit & Donahue 
2014]). There appear to be more observed systems 


with relatively low Ln a (^ 10 40 erg s 1 ) and high 
min(t coo i/tff) ( rsj 10 — 30) than our simulation for a num¬ 
ber of reasons. First, min(t coo i/tff) is measured with nar¬ 
rower radius bins in our simulations than most observa¬ 
tions, which tends to give a smaller value. Second, we do 
not include mergers in our simulations, which could trig¬ 
ger precipitation when min(£ coo i/£ff) is higher (^ 20—30). 
Third, we only simulate one system based on Perseus, a 
relatively massive cluster. It is possible that the area 
with low I/Ha (^ 10 40 erg s -1 ) will be more populated 
if we add more simulation data with systems of various 
different sizes. Lastly, the star formation rate is likely 
too high in our simulations, which is a common issue 
in numerical simulations that only include thermal feed- 

a- 

idh 


back from SNe (e.g. 
feedback model with 


Katz 1992|) 

emetic feec 


A more realistic SNe 


_ 'ack (e.g., Dalla Vecchia 

& Schaye 2008; Simpson et al. 2014| ) or including other 


physical processes such as cosmic rays (iSalem fc Bryan 


2014) and magnetic fields (Van Loo et al. 2014) would 


The peak of the Ha luminosity (Lu a > 10 4 ° erg s“ 
in our simulation appears to exceed all the observed data 
in Figure [lOj partly due to the reasons discussed above, 
and partly due to the incompleteness of the Ha sample. 
If we compare the SFR in the simulation with the SFR 
observed in the infrared (IR), we find that the two cover 
roughly the same range (Figure [TT]). The range of the 
SFR in our idealized simulation is also comparable to 
tha t fo und at l ate t imes in a cosmological simulation in 


likely suppress star formation and lower the SF'Fi. 


Dubois et al. (2010). 

Figure [lO] also shows that there is no linear anti- 
correlationbetween L^ a and the minimum ratio, as one 
might naively expect. This is because star formation 
rate is correlated with the amount of cold gas, but not 
the rate at which gas is cooling (which is correlated with 
^cooi/te ratio). Within each cycle, there is a significant 
delay up to 1-2 Gyr between when the ICM cooling rate is 
suppressed (characterized as the increase of the t coo \/t^ 
ratio) due to AGN heating and when SFR starts to de¬ 
cline at the end of the cycle due to the consumption of 
the cold gas. As a result of this hysteresis, there is no 
linear relation between SFR and the cooling rate (or the 
^cooi/te ratio), but they are related in that star forma¬ 
tion only occurs when the t coo \/t^ ratio is below some 
critical value. 

Figure [ 12 ] shows the relationship between SFR and 
the total amount of cold dense gas (with T < 10 5 K 
and number density above 10 cm -3 ). A tight correlation 
spanning a few orders of magnitude is maintained over 
the course of the simulation. The slope is generally con¬ 
sistent with what has been foun d between the SFR and 
the inferred H 2 mass for BCGs (O’Dea et al.j2008), and 
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Fig. 11.— The normalized distribution of the SFR in our stand¬ 
ard run (blue) compared with the IR observation s (red) of the 
BCGs in the ACCEPT sample ( jHoffer et al.| |2012|). For clarity, 
we assign a single value of 10 -2 M©/yr to all the upper limits in 
the observations and the times with no SFR in our simulations (the 
leftmost bar in the distribution). 



Fig. 12.— The correlation between SFR and the total amount 
of cold dense gas (T < 10 5 AT and number density above 10 cm -3 ) 
in the system. The data is sampled every 10 Myr and the color 
scheme is the same as Figure [2] 


for star forming galaxies when c omparing SFR and th e 
amount of dense molecular gas (Gao & Solomon 2004). 
A more detailed comparison with the observations can 
further constrain the star formation efficiency csf used 
in our model. We leave this to future work. 

Unlike SFR, observations have failed to find a clear 
correlation between the radio power of the A GN and the 


cooling properties of the ICM (Mittal et al. 


amount of cold gas in the system (Salome & Combes 


or the 


2003). It is possible that the radio power does not cor¬ 
relate with the a ctual mechanical power of the AGN, but 


_ P^_ 

even if it does (Cavagnolo et al. 2010), our simulation 

aFI 


shows that a clear correlation would be surprising simply 
due to the large fluctuations in the instantaneous AGN 
power: as Figure [2] and [3] show, the history of the SFR is 
rather smooth within each cycle, but the AGN power var¬ 
ies on very short timescales (< 10 Myr, the sampling time 
interval), because the accretion region is very small and 
the flow of the cold gas is not smooth. The left panel of 
Figure [T3| shows the SFR plotted against the instantan¬ 
eous Msmbh, and no correlation can be observed. How¬ 
ever, the average AGN power does correlate with ICM 


cooling and star formation as discussed in Section |3.2 
The right panel of Figure [13] shows the relation between 
SFR and the average Msmbh (calculated by smoothing 
Msmbh with a 200 Myr moving window, corresponding 
to the black line in the third panel of Figure [2]). The 
solid black line is the best linear fit: 


SFR = 13.4(Msmbh) — 20.6 


(14) 


with both SFR and (Msmbh) measured in M 0 /yr. The 
root mean square (RMS) error is 55. When averaged 
with a shorter moving window of 50 Myr, the RMS in¬ 
creases to 70. 

We also over-plot the best-fit line on the left panel for 
guidance, where the data shows large scatter and appears 
to be distributed widely, except for the lower right corner 
of the panel. This suggests that although Msmbh varies 
a lot with time, objects with very high Msmbh and very 
low SFR should be very rare. 

Significant AGN variability on short timescales has 
been used to explain the lack of a strong correlation 
between AGN activity and star formation in gala xies , as 
one would expect from the M — a relation (e.g., Hickox 


et al. 


2014 


Thacker et al. 


firms that even though the 


2014) 


star 


Our simulation con- 
brmation rate does not 


appear to correlate with the instantaneous Msmbh, a 
strong relationship does exist between SFR and the av¬ 
erage Msmbh- If we assume that the SMBH growth rate 
is ^ 1% of Msmbh, then the ratio of the mass added to 
the SMBH to that added to stars is approxim ately 10 -3 , 
consistent with the observed M — a relation (Merritt & 
Ferrarese]|2001 ). 


4.2. Thermal Conduction and Other Physics 

One important piece of physics not included in our 
simulations is thermal conduction. Conduction can be a 


very efficient heating mechan ism (|Rosner & Tucker|1989 
Narayan & Medvedev||2001|), but it is generally iounc 


Narayan & Medvedev 2001), but it is generally founc. 
that condu ction fails to stem cooling in the innermost 
~ 20 kpc ([R uszkowski & Begelman 20021 Voigt et al. 
2002; Smith et al.||2013|), or it may be sufficient to bal- 
ance cooli ng in some systems but not others with denser 
cores (e.g . Zakamska & Naray an|2003 Parrish et al. 2009 


Ruszkowski V Uh'2010 [ 201l| ). We plan to explore the 


effect of thermal conduction numerically in the context 
of AGN feedback in future work. Here we only qualitat¬ 
ively study its importance and try to address the ques¬ 
tion: does AGN feedback ever boost the core entropy to 
such a level that conduction can take over and turn a 
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Fig. 13. — Left: SFR vs. instantaneous SMBH accretion rate Msmbh (dots). No clear correlation is seen due to the large scatter. Right: 
SFR vs. < Msmbh >, the SMBH accretion rate averaged with a moving window of 200 Myr (squares). The color scheme is the same 
as previous figures (from red-early to purple-late). In both panels the solid black line shows the best linear fit to the data on the right: 
SFR = 13.4x < Msmbh > —20.6. 


cool-core cluster into a non cool-core cluster (|Guo & Oh 
2009)? 


Assuming the ICM is in hydrostatic equilibrium with 
a given gravitational potential, one can derive a family 
of quasi-steady solutions where radiative cooling is per- 
fectly balanced by thermal conduction at every radius 
(|Guo et al.||2008| |Voit||2011|). If a cluster has an entropy 
profile above this solution, conduction will make the core 
isothermal; if it is below the solution, conduction will fail 
to balance cooling and a cooling fl ow wo uld develop in the 
abs ence of other heating sources (Bertschinger & Meiksin 
1986). The dashed line in Figure |9| shows a cond uctively 
balanced solution for a 30% Spitzer conductivity (Spitzer 
1956) with the inner boundary conditions: 


n e ( 0) = 0.12 cm -3 ; T(0) = T(l) = 4keV . (15) 


Our cluster spends most of the time below the conduct¬ 
ively balanced solution. We caution that this solution 
does depend on the choice of the boundary conditions, 
and varying (n e ( 0 ), T(0)) within a reasonable range can 
shift the solution up and down by 50%. Still, given how 
briefly the cluster profile approaches or exceeds the con¬ 
ductively balanced line, it is unlikely that thermal con¬ 
duction would have enough time to transform the cluster 
core. Therefore, our simulation suggests that even with 
the help of conduction, AGN feedback seems unable to 
boost a core-core cluster to a non cool-core cluster. 

However, thermal conduction may reduce the amount 
of cold gas, and thus lower the star formation rate. On 
the other hand, the transition layer surrounding cold 
clumps may become wider due to conduction, resulting 
in an increase in the H<a luminosity. In addition, thermal 
conduction may help distribute the heat from the AGN 


jets more evenly, and it is almost certainly important in 
heating the outer parts of the cluster core 30 — 300 
kpc). Without conduction, by the end of our simula¬ 
tion, there is noticeable cooling around 100 — 300 kpc 
(Figure [5|: the temperature has decreased and the dens¬ 
ity has increased. This is also shown in Figure [6j the 
separation between the energy loss rate due to cooling 
within 100 kpc and 300 kpc increases with time. Even 
though cooling is roughly balanced by AGN heating glob¬ 
ally, without conduction the inner core gets slightly over¬ 
heated (the average cooling rate within 100 kpc declines) 
while the outer core is still under heated (cooling within 
300 kpc grows). This may also be why each burst of pre¬ 
cipitation/AGN feedback is stronger than the previous 
one. 

We have already discussed the possible impact of mer¬ 
gers on t he distribution of minimum t coo \/t^ ratio in Sec¬ 
tion 4.1 Merger-induced dynamical heating may have 
similar effects as conduction in heating the outer part of 
the cluster core around 100 — 300 kpc, and thus reduced 
the burden on the AGN. 

Another piece of physics missing in our simulations 
is magnetic fields. Even though with star formation 
we no longer see a massive cold disk that lasts indef¬ 
initely, our cluster still harbors a (still fairly massive) 
disk roughly half of the time in the simulation. /This_is 
higher than wh at observations seem to suggest (McDon-[ 


aid et al. 2010). Magneti c fields may suppo rt the cold 
filaments against gravity (Fabian et al.||2008) and delay 
the formation of a disk, and thus reduce the time the 
cluster spends in the disk phase. 

Magnetic fields and cosmic rays can also provide pres¬ 
sure support, particularly in the condensed gas. Fig- 
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ures[2] and [3] show that gas condensation at the beginning 
of each cycle always happens in a catastrophic way where 
the amount of cold gas surges up drastically Since the 
growth of individual cold clumps is mainly driven by the 
pressure difference (the cold clumps are typically under¬ 
pressured), including magnetic fields and cosmic rays in 
the simulation may slow down the growth and result in 
less bursty cooling and AGN feedback. In addition, both 
magne tic fields and cosmic rays could suppress star form¬ 
ation (|Salem & Bryan||2014| Van Loo et al.||2014|). We 
again leave this to future worlF 


5. CONCLUSIONS 

We have performed three-dimensional AMR simula¬ 
tions of an idealized, isolated cool-core galaxy cluster 
based on Perseus, focusing on the interplay between ICM 
cooling, AGN feedback and star formation over 6.5 Gyr. 
The momentum-driven AGN feedback is modeled with 
a pair of jets precessing at a small angle. The jets are 
powered by the accretion of cold gas in the close vicinity 
of the SMBH in the center of the cluster potential. We 
also include radiative cooling, self-gravity of the gas, star 
formation and stellar feedback in the form of thermal en¬ 
ergy. We perform parallel simulations to tease apart the 
roles individual physical processes play, and compare our 
simulations with the observations. The key results of this 
work are summarized below. 


1 . The radiative cooling of the ICM, AGN feedback 
and star formation are tightly linked in cool-core 
clusters. Momentum-driven AGN jets trigger the 
low entropy gas to condense out of the ICM. Some 
of the condensed gas is accreted onto the SMBH to 
power the AGN jets, but much of it forms stars. 
While AGN heating elevates the entropy of the 
ICM in the cluster center, lengthening the cooling 
time, star formation and stellar feedback gradually 
consume and erode the cold gas (typically over 1-2 
Gyr). When the cold gas vanishes, AGN feedback 
briefly turns off and allows the ICM to cool again, 
until gas with a short cooling time starts to precip¬ 
itate, which triggers the next AGN outburst. Even 
though AGN feedback is the primary heat source 
for the ICM, star formation and stellar feedback 
are more effective in depleting the cold gas, thus 
regulating the long-term AGN cycles. 

2. The range of star formation rate and the minimum 
^cooi/te ratio in the simulation generally agree re¬ 
markably well with the observations. Our cluster 
does occasionally show a higher star formation rate 
and a lower t coo \/ts ratio than observations in¬ 
dicate - this slight discrepancy could be related 
to many factors: mergers and cosmological infall, 
more effective stellar feedback, as well as the other 
physical processes not included in our simulations, 
such as conduction and magnetic fields. 

3. The instantaneous SMBH accretion rate in the sim¬ 
ulations shows large variations on short timescales 
(< 10 Myr). This explains why a clear correlation 
is not observed between the AGN radio power and 

6 Note that this could also be in part due to numer ical overcool¬ 
ing at the edges of cold clouds (Brighenti et al. |2015| 


star formation rate or the cooling properties of the 
ICM. Nevertheless, SMBH accretion and star form¬ 
ation are still related to each other, because both 
are directly linked to cooling. The star formation 
rate is tightly correlated with the SMBH accretion 
rate when the accretion rate is averaged over ^ 200 
Myr. 

4. Star formation is proportional to the total amount 
of cold gas, not the rate at which the ICM is 
cooling. The ICM responds to AGN heating very 
quickly, but the decline of the star formation rate 
occurs over a much longer period of time. As a res¬ 
ult, even though star formation only happens when 
the minimum t coo \/t^ ratio of the ICM is lower than 
some threshold (^ 20 in our simulations), there is 
no linear correlation between the two. This is con¬ 
sistent with the observations of precipitating cool- 
core clusters. 

5. The main results of this paper (the cyclical beha¬ 
vior of the cluster, the distribution of the minimum 
^cooi/te ratio, the general range of star formation 
rates, etc) are fairly robust when we exclude stel¬ 
lar feedback, or the self-gravity of the gas, or if star 
formation parameters are changed. Though these 
changes do have an effect on the strength and dur¬ 
ation of individual cycles, the ratio between SFR 
and the SMBH accretion rate, and the exact dis¬ 
tribution of SFR. 

6 . Due to its self-regulating nature, SMBH feedback 
can balance cooling with a wide range of feedback 
efficiencies, from 0.1% to 1% in numerical simula¬ 
tions without star formation. When star formation 
is included, however, even though a thermal bal¬ 
ance is still achieved with a low feedback efficiency 
of 0.1%, the star formation rate is unrealistically 
high and the system has too much cold gas for too 
long. Therefore, in the absence of other physics, a 
low SMBH feedback efficiency of 0.1% can be ruled 
out. 

7. AGN feedback balances ICM cooling in a dynam¬ 
ical way: during major AGN outbursts at the be¬ 
ginning of every cycle, there is more heating than 
cooling; as the cooling rate decreases and star form¬ 
ation consumes the fuel, heating yields to cooling 
until the beginning of the next cycle. Heating and 
cooling are balanced in an average sense, and a 
cool-core appearance is maintained despite the fluc¬ 
tuation. Our analysis shows that even with the 
help of thermal conduction, it is unlikely that AGN 
outbursts can make a cool-core isothermal. How¬ 
ever, conduction may heat the outer parts of the 
cool-core at radii of roughly 100 kpc, where AGN 
heating becomes less effective. 

Our numerical model including both momentum- 
driven AGN feedback and star formation not only solves 
the “cooling flow problem” but also reproduces the beha¬ 
vior of cool-core clusters that are largely consistent with 
observations. We do not claim, however, that we have 
found the best parameter set. Future work with a more 
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systematic parameter study and a more extensive com¬ 
parison with the observations can put better constraints 
on our model. A more sophisticated setup that includes 
mergers, conduction, and magnetic fields may also bring 
the model even closer to reality. 
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